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THEORETICAL BUCKLING LOADS OF BORON/ALUMINUM AND GRAPHITE/RESIN 
FIBER-COMPOSITE ANISOTROPIC PLATES 
by Christos C. Chamis 
Lewis Research Center 

SUMMARY 

Theoretical results are presented for designing with boron/aluminum composites 
when these composites are buckling critical. The composites are assumed to be rec- 
tangular plates with four simply supported edges. They are subjected to single and com- 
bined normal and shear loads in the plane of the plate. The plates are made from a 
unidirectional composite whose fiber direction is oriented at some angle to the load 
direction. 

The design data are presented in nondimensional form as buckling load against 
orientation angle for several plate aspect ratios. The results indicate that the buckling 
loads of boron/aluminum plates are independent of fiber orientation if the plate aspect 
ratio is greater than approximately 1. The buckling loads are moderately dependent on 
the orientation angle for plates with aspect ratios less than about 1. The buckling load 
is independent of aspect ratio in plates with aspect ratios greater than about 2. 

Comparison of buckling results for boron/aluminum composite plates and 
Thornel-75/epoxy composite plates indicates that the boron/aluminum composite plates 
resist buckling loads more efficiently than the Thor nel-7 5/epoxy composite on the basis 
of specific buckling strength. The results also indicate that the buckling loads of boron/ 
aluminum composite plates can be predicted using orthotropic theory if their aspect ra- 
tio is greater than about 1. 

The numerical algorithm used to solve the buckling problem and listing of the cor- 
responding computer code through which the results were obtained are included. 


INTRODUCTION 

Feasibility studies for the space shuttle indicate that the use of advanced fiber com- 
posite structural components can result in a considerable increase in payload in the 
shuttle system. Boron/aluminum and graphite/resin fiber composites are leading con- 
tenders for shuttle applications because these composites offer high stiffness-to-density 



and high strength-to-density ratios. Panels made from these materials will have to meet 
both material strength and buckling requirements. This report deals with a theoretical 
investigation of the buckling of flat rectangular panels made from boron/aluminum and 
Thornel-75/epoxy fiber composites. 

Several papers have appeared recently dealing with the buckling of anisotropic plates 
(refs. 1 to 8). However, d.esign data are not available for flat panels subjected to com- 
pressive loads and made from advanced fiber composites such as boron/aluminum and 
graphite/resin. The method described in reference 3 has proved efficient in buckling 
studies of boron/epoxy plates. This method is used herein to generate design data for 
boron/aluminum plates. Some data for Thornel-75 graphite/ epoxy resin plates are also 
generated for comparison purposes. Data were generated for aspect ratios of 1/2, 1, 

2, and 4. 

The panels considered are anisotropic and simply supported. They are subjected to 
combined in- plane (normal and shear) load (fig. 1). The material is a unidirectional 
composite with the fiber direction oriented at an arbitrary angle to the load direction 
(fig. 1). The analytical algorithm used is the assumed mode method in conjunction with 
the Galerkin method. A computer code was developed based on the Galerkin method, and 
the code was used to generate the theoretical design data presented herein. A brief de- 
scription of the analytical method is given in the report. All symbols are defined in ap- 
pendix A. The numerical algorithm used to solve the resulting eigenvalue problem is 
described in appendix B. Input data sample sheets with explanations are given in appen- 
dix C. A listing of the computer program with sample cases is given in appendix D. 


y 



Figure 1. - Fiber-composite plate geometry and loading - all four edges simply supported 
(aspect ratio = a lb). 


2 



I 


BRIEF DESCRIPTION OF UNDERLYING THEORY 


The underlying theory for buckling loads of anisotropic plates is described in refer- 
ence 3 with pertinent discussions in references 6 to 8. Briefly, this theory consists of 
expressing the potential energy of a plate in terms of displacement variables. Taking 
the variation of the potential energy function yields the field equation and the correspond- 
ing boundary conditions. The resulting system then is solved by the assumed mode 
technique in conjunction with the Galerkin method. 

The equation resulting after the variation of the energy function is 
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(The notation is defined in appendix A. ) The area integral represents the field equation, 
and the line integrals represent the boundary conditions. 

The assumed buckling mode described in reference 3 is represented by a Fourier 
double sine series. This mode satisfies the imposed boundary conditions, but it does 
not satisfy the natural boundary conditions if the material and structural axes do not co- 
incide. However, the mode is forced to satisfy the natural boundary conditions approx- 
imately through the Galerkin method as discussed in reference 3. 

Substituting the assumed mode in equation (1), applying the Galerkin method, and 
carrying out the algebra result in a set of linear equations which represent the eigen- 
value problem of the plate. This system is coupled for either a combination of shear 
and normal loads and/or noncoincident material and structural axes. 

The eigenvalue problem is solved by using the Power method, which is a highly ef- 
fective iterative numerical technique in seeking the largest eigenvalue of the system. 
The indicial equations which were used to generate this system and the Power method 
are given in appendix B in outline form. 
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BRIEF DESCRIPTION OF COMPUTER PROGRAM 


The numerical algorithm described in appendix B has been transformed into a com- 
puter code which is rather simple and can be generated from the information supplied 
in appendix B. A FORTRAN IV compiled listing of a computer program based on the 
algorithm in appendix B is given in appendix D, with sample cases and output. Input 
data sample sheets are given in appendix C. 

The inputs to the code are composite system identification, fiber volume ratio, 
orientation angle, plate aspect ratio, and flexural rigidities. The outputs are the num- 
ber of terms in the assumed mode series expansion required for convergence, the rela- 
tive error between the last two iteration cycles, the buckling load, and topo-plot data 
of the buckled shape of the plate normalized with respect to the largest deflection. 

The algorithm described in appendix B runs into difficulty when the shear buckling 
load of a plate is sought. In this particular case, the difficulty is bypassed by including 
normal loads which are a very small fraction of the shear load. Further discussion on 
why these difficulties arise is presented in reference 3. 


THEORETICAL DESIGN DATA 

The theoretical design data generated herein are based on the schematic illustrated 
in figure 1. In this figure, the type of loading condition, the plate geometry, and the 
fiber orientation are defined. The x-y coordinate reference system is referred to as 
the structural axes system. The fiber direction coordinate system which is located at 
the angle 9 from the structural axes system is referred to as the material axes sys- 
tem. The loading conditions are identified by N x , N , and N X y as is noted in the 
figure. 

The flexural stiffnesses required in calculating the buckling loads were calculated 
by using the computer code described in reference 9. Typical values of the elastic con- 
stants of the plate along its material axes are given in table I for boron/aluminum and 
Thornel- 75/epoxy composites with a fiber volume ratio of 0. 5. The flexural rigidities 
are computed as functions of the orientation angle using the data in table I. 
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TABLE I. - THEORETICAL UNIDIRECTIONAL COMPOSITE 


PROPERTIES AT FIBER VOLUME RATIO OF 0. 5 


Property 

Boron/aluminum 

Thornel-7 5/epoxy 

Longitudinal modulus, 

24. 2xl0 6 (35. OxlO 6 ) 

26. OxlO 6 (37. 8X10 6 ) 

N/cm 2 (psi) 

16. 8x10® (24. 3X10 6 ) 

6. 9X10 6 (1. 0X10 6 ) 

Transverse modulus, 

N/cm 2 (psi) 


0. 44X10 6 (0. 63X10 6 ) 

Shear modulus, 

8. OxlO 6 (11. 6xl0 6 ) 

N/cm 2 (psi) 



Poisson’s ratio 

0. 24 

0. 25 

o 

Density, g/cm 

2. 62 (0.095) 

1. 55 (0. 056) 

(lb/in. 3 ) 




Buckling of Boron/Aluminum and Thornel -75/Epoxy Panels 

Buckling loads for a single loading condition for panels made of boron/aluminum and 
Thor nel-7 5/epoxy are illustrated in figure 2, where the specific buckling stress has been 
plotted as a function of the orientation angle for an aspect ratio of 2. The schematic in 
the figure illustrates the type of load condition as well as the orientation angle. As can 
be seen, in this figure, boron/aluminum composites are more efficiently utilized than 



Figure 2. - Specific buckling stress of two fiber-composite plates - 
ail four edges simply supported. Fiber volume ratio, 0.5; aspect 
ratio a/b, 2. 
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Figure 3. - Buckling loads for two fiber-composite plates - all four edges simply 
supported. Fiber volume ratio, 0.5; aspect ratio a/b, 2. 


Thornel-75/resin composites in structures which are critical in buckling as measured by 
the specific buckling stress. 

Buckling load comparisons where the panels are loaded in the y-direction are il- 
lustrated in figure 3. In this figure, it can be seen that the boron/aluminum composite 
panel is considerably stronger in buckling than the corresponding Thornel-7 5/epoxy 
panel. In this plot, the nondimensional buckling load parameter is plotted as a function 
of the orientation angle for the fixed panel aspect ratio of 2. 

Buckling Loads for Individual Loading Conditions 

Design data for boron/aluminum panels which are subjected to compressive load in 
the x-direction are illustrated in figure 4 as a function of orientation for various aspect 
ratios. The important point to be noted from this figure is that the buckling load is in- 
dependent, or almost independent, of orientation angle in panels where the aspect ratio 
is approximately greater than 1. Another point to be noted is that the buckling load de- 
pends only moderately on the orientation angle in panels of aspect ratio less than 1. 

A cross-plot of figure 4 is illustrated in figure 5. The nondimensional load is 
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Figure 4. - Buckling loads for boron/aluminum composite plates, with all 
four edges simply supported, subjected to normal (N x ) compressive load. 
Fiber volume ratio, 0. 5. 
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Figure 5. - Buckling loads for boron/aluminum composite 
plates, with all four edges simply supported, subjected to 
normal (N x ) compressive load. Fiber volume ratio, 0.5. 
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plotted as a function of panel aspect ratio for various orientation angles. The insensi- 
tivity of the buckling load as a function of orientation angle in panels with aspect ratios 
greater than about 1 is clearly illustrated in this figure. The dashed line symbol is used 
to represent these curves to emphasize that buckling load values were computed only at 
the aspect ratios 1/2, 1, 2, and 4. 

Buckling loads for panels which are loaded in the y-direction only are given in fig- 
ure 6. As can be seen in this figure, the buckling loads are almost independent of the 
orientation angle in panels with aspect ratios of 1/2 and greater. The buckling load, on 
the other hand, is very sensitive to the aspect ratio in panels with aspect ratios of ap- 
proximately 2 or less. 

Buckling loads for a boron/aluminum composite panel loaded with shear only are il- 
lustrated in figure 7. The points to be noted in this figure are the following: 

(1) There is a mild buckling load dependence on the orientation angle for panel 
aspect ratios of less than about 1. 

(2) The buckling load is relatively independent of orientation angle for panel aspect 
ratios of greater than about 1. 

(3) The buckling load is very sensitive to the aspect ratio in panels with aspect ra- 
tios less than 2, and this dependence becomes rather insignificant for panel aspect ra- 
tios greater than 2. 




Figure 6. - Buckling loads for boron/aluminum composite 
plates, with all four edges simply supported, subjected 
to normal { Ny> compressive load. Fiber volume ratio, 0. 5. 


Figure 7 . - Buckling loads for boron/aluminum composite 
plates, with all four edges simply supported, subjected 
to shear (N xy ) load. Fiber volume ratio, 0. 5. 
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Buckling Loads for Two Equal Simultaneous Loadings 

Buckling loads, when the panel is loaded with equal loads in the x- and y- 
directions, are shown in figure 8. In this figure, the nondimensional buckling load pa- 
rameter is plotted as a function of orientation angle for various panel aspect ratios. 
The results in this figure show that the buckling load is slightly dependent on the 
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Figure 8. - Buckling loads for boron/aluminum composite 
plates, with all four edges simply supported, subjected 
to combined normal compressive {N x = N y ) loads. 

Fiber volume ratio, 0.5. y 


orientation angle for panel aspect ratios less than 1, and practically independent of 
orientation angle for aspect ratios greater than 1. The buckling load is sensitive to panel 
aspect ratio for aspect ratios less than or equal to 2. This dependence becomes insig- 
nificant for panel aspect ratios greater than 2. The curves of the buckling load as a 
function of the independent variables indicated in figure 8 parallel the curves of the 
buckling loads indicated in figures 4 and 6, for the individual loadings. 

Buckling loads for the case when the panel is loaded in the x-direction combined 
with shear are shown in figure 9. The curves of the buckling load for this loading con- 
dition are parallel to those of the individual cases (figs. 4 and 7). Buckling loads for 
the case when the panel is loaded in the y-direction combined with shear are illustrated 
in figure 10. The buckling load in this figure seems to be practically independent of 
orientation angle for the aspect ratios investigated. However, it is quite sensitive to 
the panel aspect ratio for aspect ratios less than 2. 
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Figure 9. - Buckling loaus for boron/aluminum composite 
plates, with all four edges simply supported, subjected 
to combined normal (N x ) and shear (N xy = N x ) loads. 
Fiber/volume ratio, 0.5. 


ASPECT 



ORIENTATION ANGLE, 0, DEG 

Figure 10. - Buckling loads for boron/aluminum composite 
plates, with all four edges simply supported, subjected to 
combined normal (N y ) and shear (N xy = N y ) loads. Fiber 
volume ratio, 0.5. 
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Buckling Loads for Three Equal Simultaneous Loadings 


Buckling loads for panels which are loaded in the x- and y-directions combined with 
shear are shown in figure 11. The schematic in this figure indicates the type of loadings 
and their respective ratios. The nondimensional buckling load is plotted as a function of 
orientation angle for various panel aspect ratios. 

Comparing corresponding curves from figures 8 and 11, it is seen that the addition 



< 
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Figure 11. - Buckling loads for boron/aluminum composite 
plates, with all four edges simply supported, subjected to 
combined normal (N y = N x ) and shear (N = N x ) loads. 
Fiber volume ratio, *0. 5. y 


of the shearing load decreases the buckling load of the panel only slightly. The point to 
be noted then is that a panel subjected to compressive loads in the x- and y-directions 
will resist almost an equal amount of shearing load for approximately the same buckling 
load. 


Buckling Loads for Two or Three Unequal Simultaneous Loadings 

Buckling loads for panels which are subjected to unequal loads in the x- and y- 
directions are shown in figure 12. The type of loading condition and respective loading 
magnitudes are illustrated in the sketch given in the figure. In this figure, the nondi- 
mensional buckling load parameter is plotted as a function of the orientation angle for 
various aspect ratios. The curves of the buckling load for this type of loading condition 
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Figure 12. - Buckling loads for boron/aluminum composite 
plates, with all four edges simply supported, subjected to 
combined normal (N v = (1/2)NJ loads. Fiber volume ratio 
0.5. Y 


ASPECT 



Figure 13. - Buckling loads for boron/aluminum composite 
plates, with all four edges simply supported, subjected to 
combined normal (N y = (1/2)NJ and shear tN™ = (1/2)N X ) 
loads. Fiber vol urn e y ratio, 0.5. 


parallel those of the cases with equal loading condition, as was previously discussed (see 
fig. 8). One additional point to be noted is that the buckling loads of panels with aspect 
ratios of greater than about 2 remain almost invariant as a function of aspect ratio when 
the orientation angle is greater than about 45°. 

Buckling data for panels loaded with unequal combined loading conditions are illus- 
trated in figure 13. The loading condition for the panel is illustrated in the schematic 
in the figure. The buckling load is plotted as a function of orientation angle for various 
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aspect ratios. A point to note in this figure is that, at some orientation angles, panels 
with aspect ratios greater than 2 could have greater buckling loads than panels with 
aspect ratios of 2. As can be seen, the buckling load for a panel with an aspect ratio of 
approximately 2 is lower than for the panel with an aspect ratio of 4 when the orientation 
angle is approximately less than 45° 

All the buckling data presented and discussed indicate that the buckling load is insen- 
sitive to orientation angle for panels with high aspect ratios. The buckling load is mildly 
sensitive to the orientation angle in panels with low aspect ratios. This observation 
leads to the important conclusion that the buckling loads of boron/aluminum composite 
anisotropic plates can be approximately determined by using classical orthotropic theory. 
This conclusion is indeed a useful one, since the buckling of orthotropic plates has been 
extensively treated in reference 5. 


POSSIBLE EXTRAPOLATIONS OF DESIGN DATA 

The design data presented and discussed were based on a fixed fiber volume ratio of 
0. 5. The data presented herein can be used to extrapolate buckling loads for plates made 
from composites with different fiber volume ratios. The results will be very close if the 
variation of the fiber volume ratio is within approximately ±20 percent of the 0. 5 value 
which was used in generating the design data. 

This extrapolation is justified since the buckling load is nondimensionalized with 
respect to both composite longitudinal modulus and thickness. It is well known that both 
composite modulus and thickness depend on the fiber volume ratio (ref. 10), and that this 
dependence is approximately linear in the fiber-volume-ratio range 0. 4 to 0. 6. In this 
sense, then, the extrapolation using the design data presented herein for fiber volume 
ratios within ±20 percent of 0. 5 should yield reasonable results. 

The computer code appended in appendix C can be slightly modified to compute the 
first natural frequency of anisotropic boron/aluminum composite panels. In reference 5, 
the vibration problem is discussed, and the analogy in computing the buckling load 
and the natural frequency is made. 


CONCLUSIONS 

The discussion of the theoretical design data presented leads to the following con- 
clusions: 

1. Design data for the buckling of unidirectional boron/aluminum panels with fibers 
oriented at any angle to the load direction have been generated and are reported herein. 
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2. Specific buckling stress comparisons showed that, in general, boron/aluminum 
composite panels are more efficient than high-modulus graphite/resin composite panels. 

3. The buckling load of boron/aluminum unidirectional panels is practically inde- 
pendent of fiber direction at high aspect ratio values. At these aspect ratios, the plate 
can be assumed to have its material axes coincide with its structural axes. Conse- 
quently, the classical buckling theory of orthotropic plates can be used to predict the 
buckling load. 

4. The buckling loads of boron/aluminum unidirectional panels are only moderately 
dependent on fiber direction at plate aspect ratios less than 1 . 

5. Boron/aluminum composite panels loaded by normal in-plane loads which are 
near the critical load can carry considerable shear load before they buckle. 

6. The buckling loads of boron/aluminum panels are practically independent of 
aspect ratio at aspect ratios greater than about 3. 

7. The buckling loads of panels with fiber volume ratios within approximately 
±20 percent of 0. 5 can be extrapolated from the design data presented herein by using 
the appropriate panel thickness and the appropriate composite longitudinal modulus. 

Lewis Research Center, 

National Aeronautics and Space Administration, 

Cleveland, Ohio, September 3, 1971, 

129-03. 
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APPENDIX A 


SYMBOLS 


a panel dimension, x-direction 

b panel dimension, y-direction 

D matrix of flexural rigidities, eq. (1) 

h thickness 

K plate stiffness matrix, elements given by eqs. (Bl) and (B4) 

L load matrix, elements given by eqs. (Bl) and (B4) 

M summation index limit on m 

m summation index 

N summation index limit on n; with subscripts, applied load 

N parameter in eq. (Bl) 

N cr buckling load 

n summation index 

p summation index 

q summation index 

W vector of coefficients in the displacement mode expansion, eq. (1) 

w displacement in z-direction 

x, y, z panel structural axes 

1, 2, 3 panel material axes 

6 variation 

e’,e" convergence tolerances 

9 orientation angle of material axes with respect to structural axes 

A buckling parameter 

p ratios defined by eq. (B5), also weight density 

a stress and matrix of stresses 

Subscripts: 

cr critical or buckling condition 

i row index 


column index 


j 

l undirectional composite property 

m buckling modes, x-direction, axes 

n buckling modes, y-direction 

r iteration cycle 

x, y, z directions associated with the respective structural axes 

1, 2, 3 directions associated with the respective material axes 
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APPENDIX B 


DESCRIPTION OF NUMERICAL ALGORITHM 


The numerical algorithm seeks the eigenvalue of the following matrix equation: 


[K] { W } = N[L] {W } 


(Bl) 


where [K] and [L] are (M X N) x (M X N) square matrices, {W } is a column matrix con- 
taining the W mn , and N is defined subsequently. The (i, j) elements of the [K] and [L] 
matrices are given by the following indicial expressions: 


m = 1(1) M; n = 1(1)N ^ 


p = 1(1)M; q = 1(1)N 


i = (m - 1)N + n; j = (p - 1)N + q J 


(B2) 


K _ n 

K ii "7 


— m 4 D n + m 2 n 2 (D 12 + 2Dg 3 ) + n 4 D 22 

a 3 ab b 3 


r p = m and q = n 


T 7i /b 2 a 2 

L n = - T 


r(B3) 


K ii = L ii = 0 P * m or q * n 


K ij = 


. 8^ 2 mnpq / 2m 2 + A D _ 87r 2 mnpq / 2n 2 + jL 

2, 2 2vl 2 2 ' 13 

a (q - n ) \p - m 


. 2 , 2 2 , \ 2 2 

b (p - m ) \q - n 


23 


f m + p and n + q odd 


L ij 


8mnpq 


(p 2 - m 2 )(q 2 - n 2 ) ^ 


K ij = L ij = 0 m + p or n + q even 


Y (B4) 
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where 


= 5s. 

D = 

3% 

D = 

N xy 

N ’ 


N 

p xy 

N 

uneter. 

Here 


N y , and 

N xy 


(B5) 


where N is the buckling parameter. Here N__, N , and N carry their algebraic 
— x y xy 

signs, whereas N is always taken as a positive quantity. 

The buckling loads are determined by finding the largest eigenvalue in equation (Bl) 
using the Power method. Several methods are available for finding eigenvalues of linear 
systems (ref. 3). The Power method is relatively easy to program and is applicable to 
nonsymmetric matrices which have real eigenvalues. The method yields the largest 
eigenvalue and the corresponding eigenvector. To apply the Power method, equation (Bl) 
is expressed in the following form: 


MW} = [K] - 1 [L] {W } 

where X = 1/N. The solution is obtained in an iterative fashion as follows: 


(B6) 


^{Wlr+lMKl'Ww}. 


(B7) 


where {W } is normalized relative to its largest element. Thus 


f , {W} 

(wL = ^ 


(B8) 


where 


A r = max ; 


i|( W r\ 


(i = 1 (1) M XN) 


(B9) 


The iteration process stops when 


r+1 


< e' 


r+1 


(BIO) 


-4 -6 

where e r is usually taken 10 < e’ < 10 . This procedure converges rapidly as long 

as the next largest eigenvalue is not close to X. The Power method runs into difficul- 
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ties when the shear buckling load is sought because for this case the magnitude of the 
two largest eigenvalues is the same. This difficulty is easily overcome by using very 
small values for one or both of the normal loads in combination with the shear load as 
is described in reference 3. 

The buckling load is obtained from the relation 


N 


X r + 1 


and 


N = — — 


= N 


V+l 


cr 


(Bll) 


Thus the smallest buckling load is obtained since is the largest eigenvalue of 

equation (B6). 

The buckling load of the plate then is determined by incrementing M and N in 
equation (Bl) and applying the Power method to compute N until the convergence 
criterion 


N(M + 1)(N + 1) < € „ ( Bl2 ) 

N(M X N) 

is satisfied. The parameter e M is usually chosen as an acceptable percentage of N. 
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APPENDIX C 


INPUT DATA 
Explanations 

Explanations appear in parentheses leadered to corresponding card. 


(composite Identification ) 


£ron 


( Nljr jypBER OF 4 0AI) CONDITIONS) 


ALU MI N 
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2 0 . 


|P,x. 
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UM AT 


py t i . 

o. 

0 

0.5 


vEJ 


OF TERMS IN THE 


(a)) 

(y-EDGE 




ONE SET) 


IGIDITI^S ) 
D ! t 2 
4 . 0 * 6 * 


8 9 (O II 12 13 14 15 16 17 10 


NASA-C-836 IREV 9-14-59) 


. P *y, 

o. 

1.0 
0 . 0 


D 1 3 

0 . o’ 


(CONTINUE WITH 


SERIES 


fle: 
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DIMENSION 


31 


4 . 


19 20 21 22 23 24,25 26 27 28 29 3q 


06 


37 38 39 40 61 42 


63 64 45 66 67 < 


IXURAL RIGIDITIES 


(b)) 

y* 

1.0 

i . 6 

0.5 


d 22 ( 

7 .0 1 


31 32 33 34 35 36 


k9 50 51 52 53 54 


■37 38 39 60 41 62 


23 , 


OTHER 


(NUMBER 


0* 

5 


ANGLES) 


[63 66 65 66 67 68 


55 56 57 58 59 60 


LOADINi 


61 62 63 66 65 66167 68 69 70 71 


CONDITION 


49 50 51 52 53 56 55 56 S7 S8 59 60 


(or: 


61 62 63 66 65 66 


JjENTATIOrf 

) 


:ts 


72 73 74 75 76 77 78 79 81 




67 68 69 70 71 72 




ANGLE) - 


73 74 75 76 77 78 79 8 
CD <,729 
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Sample Input Data for 0. 02-lnch-Thlck Boron/Aluminum Plate 
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>1 ^ ^ 


APPENDIX D 


COMPUTER OUTPUT 


C.... BUCKLING OF ANISOTROPIC PLATES 
CCMMCN /EIGVEC/ W(81) 

CCMMCN /ARFIL1/ D ( 3 , 3 )./ ARF IL 2/RHO( 3 ) /ARR AY/X ( 81 9 81 ) 

REAL N1,?^2,LAMB2 

DIMENSICA VdOvlOI't TflO), RST ( 20, 3 ) t T ITLE { 3 ) 

DATA P I / 3 • 141 5927/ 

MOCE= 1 

REAC ( 5 , 23 ) TITLE 
RE AC (5,12) NP, MU 
REAC (5,13) A,B 

READ (5,14) ( (RST( I»JV,J=1,3),I=1,NP) 

1 REAC (5,25) ( <D( UJ),J=1,3), 1=1, 3 ), THETA 
DC 11 I 1*1 ,NP 

WRITE (6,15) 

DC 2 J= 1 , 3 

2 RHC(J);=RST (II , J) 

N 1 = 0 • C 

M= 1 

CALL TIMLFT ( Q1 ) 

3 

CALL MFILL1 (M,A,B) 

CALL MINN/A ( MM, K S I G ) 

IF (KSIG.EC.2) GO TO 10 
CALL MFILL2 t M, A , B ) 

CALL EIGFIN ( MM, L AMB2 ) 

IF (KSIG.EC.3) I FL AG =1 
N2 = 1 *0/ L A M 8 2 
EPS=ABS( (N2-N11/N2) 

IF (EPS.LE.1.0E-4) GC TO 5 
IF (R.EC.MU) GO TO 4 
SPS= EPS 

m=m + 1 

N1 = N2 
GC TC 3 
WRITE (6,16) 

IF (IFLAG.EQ.l) WRITE (6,17) 

CALL T I ML FT ( Q2 ) 

IC=(Cl-C2)/60. 

WRITE (6,24) TITLE, THETA 
WRITE (6,18) M, S PS , EPS , 1 0 
VM A X =0 • 0 
DC 7 1=1,10 
DC 7 J=l,10 
VI I • J )= 0 * 0 
DC 6 K=1,V 
DC 6 L= 1 , M 
IJ=(K-1)«m+l 

V( I , J)=V ( I , J )+'W( I J ) *S IN ( FLOAT ( K* l )*P I/10.0)*5IN(FL0AT(L*J)*PI/10.) 
VMAX=AMAX 1 (VMAX, ABS(V( I, J) ) ) 

DC 8 1=1,10 
DC 8 J=l,10 
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I 


8 V( I t J ) = V ( I , J ) / VM A X 

WRITE (6,19) A,8,N2 f MODE,RHO,D 

DC 9 1=1,10 
J=1C-I+1 
T < I ) =0. 0 

9 WRITE (6,20) J,T < I ) , < V< K, J ) , K=l, 10) 

1=0 

XI=C.O 

WRITE (6,20) ItXIv (T( J) 9 J-lt 10) 

WRITE (6,21) 

GO TC 11 

10 WRITE (6,22) 

11 CONTINUE 
GO TC 1 

C 

12 FORMAT (212) 

13 FCRN'AT (2F10.4) 

14 FCRVAT (9F8.3) 

15 FCRVAT ( 1 H 1 ) 

16 FCRVAT (23H MAXIMUM INDEX ATTAINED) 

17 FORMAT (41H ILL-CONDITIONED MATRIX - RESULT IN DOUBT) 

18 FCRVAT C:X,2HM=, 1 1 , 5X , 20HPR EV I OUS REL. ERROR = , 1 D E 10 • 3 , 5X , 1 6HL A S T R 
1EL. ERR0R=,1PE10.3,5X^33H TIME REQUIRED FOR THIS CASE WAS ,I3,8H S 
2ECCNCS ) 

19 FORVAT (///7X,1HA,13X„1HB, 1 1 X , 3HNCR , 1 1 X , 4HM0DE , 9X , 5 HRHO- X , 9 X , 5HRH0 
1-Y ,9X ,6HRHC-XY// 3 ( 2 X , l P E 10 . 3 , 2X ) , 6X , 1 1 , 7X t 3 ( 2X , 1PE 10. 3 , 2 X ) / / /37H T 
2 HE C-ARRAY COLUMNWISE IS AS FOL L OWS / / 9 ( 2 X , 1 P E 10 . 3 , 2 X ) / / / /50 X , 1 3H BU 
3CKLEC SHAPE// 6X, 1HY) 

20 FORVAT ( / / 5 X , I2,11(5X , F 5*2)) 

21 FORVAT (//l5X,lH0*9X,lHl,9X,lH2,9X,lH3,9X t lH4,9X,lH5,9X f lH6,9X,lH7 
1 ,9X f 1H8 ,9X , 1H9,8X t 2H10, 6X, 1HX ) 

22 FCRVAT (37H MATRIX SINGULAR - PROCESS TERMINATED) 

23 FCRVAT (3A6) 

24 FCRVAT ( /20H COMPOSITE SYSTEM = ,3A6,12H AT THETA = ,F5.1) 

25 FCRVAT (10FR.3) 

ENC 


SUBROUTINE VF ILL 1 ( M V ,. A , B ) 

CCVVCN /AP.FIL1/ 0(3,3) 

CCVVCN /ARRAY/ X (81 ,81) 

INTEGER F,C 

DATA P I / 3 • 141 592 7/ 

DC 3 V=l f VV 

DC 3 N= 1 , V v 

DC 3 P=1,W 

DC 3 C = 1 , V v 

I^= ( V-l ) *NM + N 

J= ( P-1 ) *NM+C 

IF ( I .NE. J ) GO TO 1 

IF (P.NE. V.CR.Q.NE.N) GO TO 2 

X ( I , J ) = PI **4* ( B+FLOAT (M**4 ) *C( 1 , 1 ) / A**3 + 2,0*FLOAT( V*M*N*N) * (0< 1 ,2 ) 
1 +' 2 . C * C ( 3,3 ) ) / (A*P)+A*FLnAT(N**4)*D( 2, 2 ) / B * * 3 ) / 4 * 0 
GC TC 3 
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o o O \J> o -p* 


1 IF (MCD< (M+P)*(N+0) ,2).EQ.O) GO TO 2 

X < I , J J=-8 .0*PI**2*D(lf3)*FL0AT(M*N*P*Q*(2*M*M/(P*P-M*M>+l)/(Q*0-N* 
IN) ) /A **2-8 # 0*PI**2*D( 2, 3 )*FLOAT( M*N*P*Q*( 2*N*N/ (Q*Q-N*N)+1) / (P*P-M 
2*M) ) /B**2 
GC TC 3 

2 XU,J> = 0.0 

3 CONTINUE 
RETURN 
END 


SUBROUTINE MINVA ( M , KS I G ) 

C... ..MATRIX INVERSION BY GAUSS-JORDAN ELIMINATION 
COMMON / £ R R AY / X ( 8 1 1 8 1 ) 

DIMENSION Y(81,81), Z(2»8l,81)t K(81) 
EQUIVALENCE (Y,Z) 

DOUBLE PRECISION A,R t PR0D,Y 

DATA CELT , EPS t LOOPS /O, OOOI f 1.0E-8, 2/ 

KSIG=C 

IF (M.NE.l ) GO TC 3 
PRC C = X ( 1,1 ) 

IF (PROD. NE. 0.0) GC TO 2 

1 K$IG=KS IC+2 
RETURN 

2 XU, 1) = 1.0/PR0D 
RETURN 

3 PRO C= 1 • 0 

MM=v-l 

DC A 1=1, M 
K( I) = I 
DC A J= 1 1 M 
YU , J) = XII , J) 

..BEGIN BY FINDING LARGEST PIVOTAL ELEMENT 
DC 9 I = 1 , v 
A = C.O 

DC 5 J= I , M 

IF (ABS(V( J,in.LE.A) GO TO 5 
A~ARS(Y(Jt I) ) 

L= J 

CONTINUE 

IF (A. EC. 0.0) GO Tn 1 
..REARRANGE ROWS AND ORDER ARRAY 
N-K (L) 

K { L ) = K ( I ) 

K(I )=N 
DC 6 J- 1 » M 
A=V ( I *J ) 

YU,Jt=Y(L,J) 

Y ( L , J) = A 

..RECLCE PIVCTAL ROW 
A = Y ( I » 1 ) 

DC 7 J= 1 t M M 
7 Y( I , J) = Y(I , J+l) / A 

Y ( I,M)=1.Q/A 
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I II I 


C REDUCE REMAINING ROWS 

DC 9 L= 1 » M 
IF (L.EC. I ) GO TO 9 
A-YIU1 ) 

DC 8 N=1,MM 

Y(L,N»=Y(L,N+i)-A*Y( I,N) 

8 IF (AeS(Y(LtN)).LT.(ABS(Y(L,N+l) )*EPS») Y(L,N)=O.C 
Yt L , M ) = -A *Y ( I *M) 

9 CONTINUE 

C UNSCRAMBLE INVERTED MATRIX 

DC 13 I = 1 , M 



IF 

f K ( I ).EC.I ). 

GC 

TO 

13 


DC 

10 J=I f P 




10 

IF 

( K ( J ) • EC • I ) 

GC 

TO 

11 


GC 

TC 1 




11 

DC 

12 L = 1 t Y 





A=Y(L,I ) 





Y ( L 

,I)=Y(L,J) 




12 

YCL 

<1 

II 

-3 





K( J)=K( I ) 




13 

CONTINUE 





C OBTAIN ERROR MATRIX 

DO 18 N=l, LOOPS 
TES T = G. 0 
DC 15 I = 1 » M 
DC 15 J = 1,M 
R = 0 .0 

DC 1 A L = 1 , M 

14 R=R-7(1,L,J)*X(I ,L) 

IF (I.EO.J) R = R+ 1.0 
TEST = AMAX II TEST , ABS(R).) 

15 Z(2,UJ)=R 
DC 17 1 = 1 ,M 
DC 17 J= 1 » M 
A = 0 • 0 

DC 16 L = 1,M 

16 A=A+Z(1,I,L)*Z(2,L,J) 

17 Z ( 1 , I ,J )=Z ( 1, I, J )+A 

IF ( TEST.LE.DELT ) GO TO 19 

18 ctntinue 

KS IC=KS 10+3 

C.... .TRANSFER FINAL INVERSE 

19 DC 2C 1=1 ,M 
DC 2C J = 1 » M 

20 Y< I,J) = Z«1, I,J) 

RETURN 

ENC 


SLRPCLTINE MF ILL 2 (MM t A,B) 
COMMON / A R F I L 2/ RH0(3) 

COMMON /ARRAY/ X (81,81) 

DIMENSION 0(81,81), 0(81,81) 

REAL L 

INTEGER F,C 

DATA PI /I . 1415927/ 
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-4 O' 


min i 


DC 3 P=1,PP 

DC 3 N=1,PP 

DC 3 P=1,PP 

DC 3 C=1,PP 

1= (P-1) *PI«+N 

J= (P-1) *M+C 

IT ( I.NE. J) GO TO 1 

IF (P.NE.M.OP.Q.NE.N) GO TO 2 

L ( I » J > = -F I **2*( B*FLOAT( Pl*PI)*RHO( 1 ) / A+ A*FLO AT ( N*N ) *RHO ( 2 ) /B ) /4. 0 
GC TC 3 

1 IF <PCD( (P» + P)*(N+Q) ,2).EQ.OI GO TO 2 

H I » J I=-RHC (3 )*FLOAT( 8*M*N*P*Q/ < ( P*P-M*M ) * ( Q*Q-N*N ) ) ) 

GO TC 3 

2 Ltl,J>=0.0 

3 CONTINUE 
NN=PP*MP 

DC 5 1=1 ,NN 
DC 5 J=1,NN 
SUP = C .0 
DC 4 I J = 1 » NN 

4 SUP=SUM+X ( It I J)*L( IJ, J) 

5 C(I,J).= SLP 
DC 6 1=1, 'NN 
DC 6 J= 1 » A N 

6 X( I , J )=C ( I , J ) 

RETURN 

END 


SUBROUTINE EIGFIN { PP..L AMB21 
REAL LAPP 1 , L AM82 
CCPPCN /EIGVEC/ W t B 1 ) 

CCPPCN /ARRAY/ X(81,8l> 

DIPENSION V < 8 1 ) 

LAVP1=0.C 
W( 1 )=0.5 

IF (PP.EC.l) GO TO 2 
DC 1 I = 2, PR 

1 W ( I ) =C . 5 

2 DC 4 1 = 1, PR 
SUP=C.D 

DC 3 J=1,PP 

3 SUP=SUM+X(I,J)*W(J) 

4 VU) = SUP 

L AN B2 = ABS ( V( 1 ) ) 

IF (PP.EC.l) GO TO 6 
DC 5 1=2, PP 

5 LAPe2 = AVAXl(LAMB2, ARS(V( I) ) ) 

DC 7 1 = 1, PP 
MI ) = V ( I ) / L AMR2 

IF (APS ( (LAPR 2-L APR1I/LAMB2) .LE.1.0E-4) RETURN 
L ANp ItL AN B 2 
GC TC 2 
END 
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MAXIMA IK'CEX ATTAINED 


C0MPC3ITE SYSTEM = BORON /ALUMINUM AT 

P=;5 PREVIOUS REL. ERROR* 1.554E- 

THETA = 30.0 

04 LAST REL. 

ERROR = 3.944E 

-04 

TIME REQUIRED FOR 

THIS CASE WAS 

11 SECONDS 

A 


B 

NCR 

MODE 

RH0-X 

RH0-Y 


RH0-XY 



. 2.CC0E+01 

1 

.OOOE+OI 

8.111E+00 

l 

-1 • C0OE+OO 

0. 

0 

• 



THE C- ARRAY 

COLUMNWISE IS AS 

FOLLOWS 








2.219E+C1 

4 

•490E+00 

1. 870E+00 

4. 490E+00 

1* 845E+01 

1 • 370E+00 1 

• 870E+00 

1.370E+00 

8.200E+00 

V 




BUCKLED 

SHAPE 






1C 

0. 

-0.00 

-0.00 

i 

o 

. 

o 

0 

1 

o 

. 

o 

o 

0.00 

0.00 

0.00 

0.00 

0.00 

-0.00 

9 

0, 

0.20 

0.31 

0.26 0.12 

-0.04 

-0.14 

-0.18 

-0.15 

-0.09 

0.00 

8 

0. 

0.39 

0.58 

0.50 0.22 

-0.08 

-0.27 

-0.34 

-0.28 

-0.16 

0.00 

7 

0. 

0.55 

0.80 

0.67 0.29 

-0.11 

-0.38 

-0.46 

-0.39 

-0.22 

0.00 

6 

0. 

0.65 

0.95 

0.79 0.33 

-0.14 

-0.45 

-0.54 

-0.45 

-0.25 

0.00 

5 

0. 

0.69 

1.00 

0.82 0.34 

-0.16 

-0.47 

-0.56 

-0.47 

-0.26 

0.00 

4 

0. 

0.66 

0.95 

0.77 0.31 

-0.16 

-0.45 

-0.54 

-0.44 

-0.24 

0.00 

3 

0. 

0.57 

0.81 

0.65 0.25 

-0.14 

-0.39 

-0.45 

-0.37 

-0.20 

0.00 

2 

0. 

0.42 

0.59 

0.46 0.18 

-0.11 

-0.28 

-0.33 

-0.27 

-0.14 

0.00 

1 

0. 

0.22 

0.31 

0.24 0.09 

-0.06 

-0.15 

-0.17 

-0.14 

-0.07 

0.00 

0 

0. 

0. 

0. 

0. 0. 

0. 

0. 

0. 

0. 

0. 

0. 


c 

1 

2 

3 4 

5 

6 

7 

8 

9 

io x 



M COMPCSITE SYSTE M = BORON/ ALUM I NUM AT THETA = 30.0 


M3 

PREVIOUS REL • ERRORS 8.298E- 

04 LAST 

REL. 

ERROR= 3.563F- 

06 TIME REQUIRED FOR 

THIS CASE WAS 

5 SECONDS 

A 


B 

NCR 

MODE 


RHO-X 

RH0-Y 


RHO-XY 



2.CC0E+C1 


1.000E+01 

2 • 9R6E+00 

1 


0. 

-1 • OOOE+OO 

0 

• 



the c-array 

COLUMNWISE IS AS 

FOLLOWS 









2.219E+CI 


4.490E+00 

1.870E+00 

4.490E+00 

1 • 845E+01 

1 • 370E+00 

1 

• 870E+00 

1.370E+00 

8.200Ef00 





BUCKLED 

SHAPE 






Y 












10 

0. 

-0.00 

-0.00 

-0.00 

-0.00 

o 

o 

. 

0 

1 

-0.00 

-0.00 

-0.00 

-0.00 

0.00 

9 

0. 

0.09 

0.18 

0.24 

0.29 

0.31 

0.30 

0.25 

0.18 

0.10 

o 

o 

. 

0 

1 

e 

0. 

0.18 

0.34 

0.47 

0.55 

0.59 

0.56 

0.48 

0.35 

0.18 

-0.00 

7 

0. 

0.25 

0.47 

0.65 

0.76 

0.81 

0.77 

0.66 

0.40 

0.25 

-0.00 

6 

8. 

0.29 

0.55 

0.76 

0.90 

0.95 

0.90 

0.77 

0.56 

0.29 

-0.00 

5 

0. 

0.31 

0.59 

0.81 

0.95 

1.00 

0.95 

0.81 

0.58 

0.31 

-0.00 

4 

0 . 

0.30 

0.57 

0.78 

0.91 

0.95 

0.90 

0.76 

0.55 

0.29 

o 

o 

. 

0 

1 

3 

0 . 

0.26 

0.49 

0.66 

0.78 

0.81 

0.77 

0.65 

0.47 

0.24 

-0.00 

2 

0 . 

0.19 

0.36 

0.49 

0.57 

0.59 

0.55 

0.47 

0.34 

0.18 

o 

© 

. 

o 

t 

1 

0 . 

0.10 

0.19 

0.26 

0.30 

0.31 

0.29 

0.25 

0.18 

0.09 

-0.00 

C 

0 . 

0 . 

0 . 

0 . 

0 . 

0 . 

0 . 

0. 

0 . 

0. 

0. 


c 

1 

2 

3 

4 

5 

6 

7 

8 

9 

10 X 



COMPOSITE SYSTEM = BORON/ ALUM I NUM AT 

M*3 PREVIOUS REL. ERROR= 8.298E* 

THETA = 30.0 

04 LAST REL. 

ERROR 3 5.726E 

*06 TIME REQUIRED FOR 

THIS CASE WAS 

6 SECOND S 

A 


B 

NCR 

MODE 

RHO-X 

RH0-Y 


RH0-XY 



2.CC0E+01 

1 

•OOOE+Ol 

2. 389E + 00 

1 

- l • OOOE+OO 

-l . OOOE+OO 

0 

• 



THE C-ARRAY 

COLUMNWISE IS AS 

FOLLOWS 








2.219E+01 

4 

.490E+00 

1 . 870E + 00 

A. A90E+00 

1 . 8A5E+01 

1 • 370E+00 

1 

• 870 E+00 

1 .373E + 00 

8.200E+0C 

Y 




BUCKLED 

SHAPE 






1C 

0 . 

-0.00 

-0.00 

o 

o 

0 

1 

o 

0 

o 

1 

-0.00 

-0.00 

-0.00 

-0.00 

-0.00 

0.00 

9 

0 . 

0.10 

0.19 

0.26 0.30 

0.31 

0.29 

0. 24 

0.17 

0.09 

-0.00 

8 

0 . 

0.19 

0.36 

0.A9 0.57 

0.59 

0.55 

0.46 

0.33 

0.17 

-0.00 

7 

0 . 

0.27 

0.50 

0.68 0.78 

0.81 

0.75 

0.63 

0.45 

0.24 

-0.00 

6 

0 . 

0.31 

0.59 

0.80 0.93 

0.95 

0.88 

0.74 

0.53 

0.27 

-0.00 

5 

0 . 

0.33 

0.63 

0.85 0.98 

1.00 

0.93 

0.77 

0.55 

0.29 

“0.00 

4 

0 . 

0.32 

0.60 

0.81 0.93 

0.95 

0.88 

0.73 

0.52 

0.27 

-0.00 

3 

0 . 

0.28 

0.52 

0.69 0.79 

0.81 

0.7A 

0.62 

0.44 

0.23 

“0.00 

2 

0 . 

0.20 

0.38 

0.51 0.58 

0.59 

0.5A 

0.45 

0.32 

0.16 

-0.00 

1 

0 . 

0.11 

0.20 

0.27 0.3C 

0.31 

0.28 

0.23 

0.17 

0.09 

-0.00 

C 

0 . 

0 . 

0 . 

0 . 0 . 

0 . 

0 . 

0 . 

0 . 

0 . 

0 . 



NAXIMM INCEX ATTAINED 


CCMPC SITE SYSTEM = BOR ON / ALUM I NUM AT 

M=i5 PREVIOUS REL. ERROR= 1.910E- 

THETA 

02 

= 30.0 

LAST REL. 

ERROR= 1.710E 

-03 

TIME REQUIRED FDR 

THIS CASE WAS 

ll SECONOS 

A 


B 

NCR 

MODE 

RHO-X 

RHO-Y 


RHD-XY 



?• C CCE + 0 1 


1.000E+01 

6.185F+00 


1 

-1 • OOOF+OO 

0. 

-1 

• OOOE+OQ 



THE C-ARRAy 

COLUMN’ WISE IS AS 

FOLLOWS 









2.219E+01 


A.A90E+00 

1.87CE+00 

A.A90E+00 

1 • 8A5E+01 

1 • 370E+00 1 

• 870E+00 

1 • 370E + 00 

8.200E4-00 

Y 





BUCKLED 

SHAPE 






1C 

0. 

-0.00 

-0.00 

-0.00 

-0.00 

-0.00 

0.00 

0.00 

0.00 

0.00 

0 
o 

• 

D 

1 

9 

0. 

0.10 

0.22 

0.32 

0.30 

0. 16 

-0.06 

-0.2 A 

-0.31 

-0.22 

0.00 

8 

0. 

0.21 

0.A5 

0.61 

0.55 

0.25 

-0.17 

-0.50 

-0.58 

-0.39 

0.00 

7 

0. 

0 • 3 A 

0.66 

0 • 8 A 

0.70 

0.23 

-0.3A 

-0.73 

-0.79 

-0.50 

0.00 

6 

0. 

0.A5 

0.83 

0.97 

0.72 

0.12 

-0.53 

-0.91 

-0.89 

-0.53 

0.00 

5 

0. 

0.53 

0.92 

0.99 

0.63 

-0.03 

-0.68 

-1.00 

-0.90 

-0.51 

0.00 

A 

0. 

0.55 

0.91 

0.90 

0.A7 

-0.19 

-0.76 

-0.97 

-0.81 

i 

o 

• 

u> 

0.00 

3 

0. 

0.51 

0.80 

0.71 

0.29 

-0.29 

-0.72 

-0.R3 

-0.6A 

-0.32 

0.00 

2 

0. 

0 . AO 

0.59 

0.A7 

0.13 

-0.29 

-0.57 

-0.60 
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